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Abstract —Reachability analysis is important for studying 
optimal control problems and differential games, which are 
powerful theoretical tools for analyzing and modeling many 
practical problems in robotics, aircraft control, among other 
application areas. In reachability analysis, one is interested in 
computing the reachable set, defined as the set of states from 
which there exists a control, despite the worst disturbance, that 
can drive the system into a set of target states. The target 
states can be used to model either unsafe or desirable configu¬ 
rations, depending on the application. Many Hamilton-Jacobi 
formulations allow the computation of reachable sets; however, 
due to the exponential complexity scaling in computation time 
and space, problems Involving approximately 5 dimensions 
become intractable. A number of methods that compute an 
approximate solution exist in the literature, but these methods 
trade off complexity for optimality. In this paper, we eliminate 
complexity-optimality trade-offs for time-invariant decoupled 
systems using a decoupled Hamilton-Jacobi formulation that 
enables the exact reconstruction of high dimensional solutions 
via low dimensional solutions of the decoupled subsystems. Our 
formulation is compatible with existing numerical tools, and we 
show the accuracy, computation benefits, and an application of 
our novel approach using two numerical examples. 

I. Introduction 

Optimal control problems and differential games have 
been extensively studied [1], [2], [3], [4], and have received 
growing interest in the recent past. These powerful theoreti¬ 
cal tools allow us to analyze a variety of real world problems, 
including path planning, collision avoidance, safety verifica¬ 
tion, among other applications in robotics, aircraft control, 
security, and other domains [5], [6], [7], [8]. 

In an optimal control problem, one aims to drive a con¬ 
trolled dynamical system into a set of states called the target 
set; depending on the application, the target set can model 
the set of either desirable or undesirable configurations. In a 
reachability framework, one aims to determine the backwards 
reachable set, defined as the set of states from which a control 
exists to drive the system into the target set. Differential 
games involve two adversarial players (Player 1 and Player 
2). Player 2 seeks to drive a system to a target set, while 
Player 1 seeks to prevent Player 2 from doing so. One again 
aims to determine the backwards reachable set, which in this 
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case is defined as the set of states from which a control from 
Player 2 exists to drive the system into the target set, despite 
the optimal adversarial control from Player 1. 

Reachability is an effective way to analyze optimal con¬ 
trol problems and differential games because it provides 
guarantees on system performance and safety. Reachability 
problems involving one player can be posed as a minimum 
(maximum) cost game where the player minimizes the mini¬ 
mum value over time of some cost function representing the 
proximity to the target set. In the case of a differential game. 
Player 1 maximizes the minimum cost over time, while 
Player 2 minimizes it. [1] has shown that the backwards 
reachable set can be obtained by solving a Hamilton-Jacobi 
Partial Differential Equation (HJ PDE) with a terminal con¬ 
dition specifying the target set. Many similar formulations of 
the backwards reachability problems also exist [9], [10], [11]. 
HJ reachability has been successfully used to solve problems 
such as aircraft collision avoidance [1], automated in-flight 
refueling [12], and reach-avoid games [13], [14]. 

The techniques for computing backwards reachable sets 
via solving an HJ PDE are very flexible and can be applied 
to a large variety of system dynamics when the problem 
dimensionality is low. Furthermore, many numerical tools 
have been developed to solve these equations, making the 
HJ approach practically appealing [15], [16], [17], [18]. For 
higher dimensional problems, various techniques such as 
those involving projections [19], [20], approximate dynamic 
programming [21], and occupation measure theory [22] 
have been proposed. While these approximation techniques 
alleviate the computation complexity, they give up optimality 
and sometimes give overly conservative results. 

This paper resolves the complexity-optimality trade-off for 
time-invariant systems with decoupled dynamics. We present 
a decoupled formulation of HJ reachability for decoupled 
systems, defined in (|^. By considering the decoupled com¬ 
ponent separately and solving lower dimensional HJ PDEs 
for each subsystem, we reduce the computation time and 
space complexity substantially. Our approach also exactly 
recovers the solution to the original, high dimensional PDE. 

II. Problem Formulation 

Consider a differential game between two players de¬ 
scribed by the time-invariant system 

z = f{z, u, d), almost every t G [—T, 0], (1) 

where z G M" is the system state, u G Id is the control of 
Player 1, and d G V is the control of Player 2. We assume 
/ : K” X U X V ^ M” is uniformly continuous, bounded. 


and Lipschitz continuous in z for fixed u, d, and the control 
functions u{-) G U,d(-) G D are drawn from the set of 
measurable function^ As in [1], [23], [24], we allow Player 
2 to only use nonanticipative strategies 7, defined by 


to compute the backwards reachable set, V(t), in the low¬ 
dimensional space K"* of each of the decoupled components 
Xi as opposed to in the full system state space K". V{t) is 
defined as 


7 e r := {A/" : U -)> D I 

u{r) = u{r) for almost every r G [f, s] 

^ A/"[u](r) = A/’[M](r) for almost every r G [t, s]} 

( 2 ) 

We further assume that the system is a decoupled system. 

Definition 1: Decoupled system. A system Q is a decou¬ 
pled system if it can be split into N components, denoted 
where z = (xi,... ,Xm), that satisfy the following: 


V{t) := {z S K" I 37 G r such that 

Vu(-) e U,3s G [f,0],C/(s;z,f,u(-),7M(-) G -C)} 

(7) 

Remark 1: One may have noticed that if £ = [J^ Ci, one 
would be able to simply find Vi (t) by solving ([^ with Ci as 
the target set, and then obtain V{t) = [JiVi(t). However, it 
is crucial to observe that we are interested in the case where 
£ = Ci, in which a simple union of Vi{t) would not yield 
the correct reachable set V(t). 


Xi = fi{xi,Ui,di), almost every t G [-T,0], 


where Xi G M"’ is ith component of the full state, Ui G Ui 
is ith component of the control of Player 1, and di G Di 
is *th component of the control of Player 2. Based on this 
assumption and the assumptions on /(•,•,•), u(-), d(-), we 
have that fi : xUiXVi -G M"’ is uniformly continuous, 

bounded, and Lipschitz continuous in Xi for fixed Ui,di and 
Ui{-),di{-) are measurable. Note that = n. 

Denote system trajectories, which are solutions to ([T]), as 

^/(s; z, t, u{-), d{-)) : [t, 0] M”. (4) 


satisfies initial conditions ^f{t;z,t,u{-),d{-)) = z and 
the following differential equation almost everywhere 


—^fis;z,t,u{-),d{-)) = f{^f{s;z,t,u{-),d{-)),u{s),d{s)) 

(5) 

In our differential game, the goal of Player 2 is to drive 
the system into some target set £, and the goal of Player 1 
is to drive the system away from it. The set £ is represented 
as the zero sublevel set of a bounded, Lipschitz continuous 
function Z : K" —)• K, £ = {z G K" | l{z) < 0}. 

Such a function always exists, since we can choose l{-) 
to be a signed distance function; we call l{-) the implicit 
surface function representing the set £. In accordance with 
our decoupled dynamics, we assume that I can be represented 
as a maximum of N bounded, Lipschitz continuous functions 
li : —>■ M, l{z) = l{xi,... ,xn) = maxpi(xi) where 

li{xi) are implicit surface functions representing Ci so that 
z G C ^ Xi G Ci Mi. Note that with the definition of l{z) 
and li{xi), we have that 


£=n^- (6) 

i 

Given the decoupled system and the target set £ in 
the form (|^ represented by Z( ), our goal in this paper is 

*A function f ■. X ^ Y between two measurable spaces (X, Sx) and 
(y, Sy) is said to be measurable if the preimage of a measurable set in 
y is a measurable set in X, that is: VP S S Sx, with 

Sx, Sy (T-algebras on X,Y. 


III. Solution 

A. HJ Reachability: Full Formulation 

In [1], the authors showed that the backwards reachable set 
V(f) can be obtained as the zero sublevel set of the viscosity 
solution [25] V {t, z) of the following terminal value HJ PDE: 

DtV{t,z) -I-min{0, maxmin£) 2 L(£ z) • f{z,u, d)]} = 0 

u^lA d^T> 

V{0,z) = liz) 

( 8 ) 

from which we obtain V{t) = {z G K" | V(t,z) < 0} 
from the bounded, Lipschitz function V{t, z) that is also 
continuous in both z and t [23]. 

[1] and similar approaches, such as [9], [10], [11], are 
compatible with well-established numerical methods [15], 
[16], [17], [18]. However, these approaches become in¬ 
tractable quickly as the dimensionality of the problem n 
increases. Numerically, the solution V{t,z) is computed on 
a grid, and the number of grid points increases exponentially 
with the number of dimensions. 

Decoupled dynamics allow for tractable or faster compu¬ 
tation of reachable sets in the individual decoupled com¬ 
ponents. Some authors [19], [20] have proposed methods 
for combining or stitching together these reachable set com¬ 
ponents into the full reachable set. These methods work 
reasonably well, but introduce conservatism in various ways. 
In the next subsections, we will provide a method for 
combining solutions to the lower dimensional HJ PDEs to 
construct the exact full solution in the original HJ PDE. 

B. HJ Reachability: Decoupled Formulation 

Observe that can be viewed as an equation involving 
two cases. Depending on which of the arguments in the outer¬ 
most minimum is active, i) becomes one of (|^ or 

DtV{t,z)=0 (9) 


DtV{t, z) + maxmin£>zL(f, z) ■ f{z,u,d) = 0 (10) 


This motivates us to define the following sets 
which characterize which of the outer-most minimum oper¬ 
ation is active in I®). 


Theorem 1: The solution to (|^ for a decoupled system 
with dynamics ([^ and terminal condition l{z) = max^ li{xi) 
is given by 


= {z S M" I maxminZ? 2 l/(f, z) • f{z,u,d) > 0} 

u^U d^V 

d- 2 {t) = {z e M" I max min DzV{t,z) ■ f{z,u,d) < 0} 
ueu d^v 

( 11 ) 

Note that Fi{t) is the complement of F 2 {t), d-i{t) = 
for all time. At a given t, in V(t,z) sat¬ 
isfies in J^ 2 {t), V{t,z) satisfies We now show 

an important property of and J- 2 {t) in the Lemma 

and Corollary below. These will be used to show that our 
proposed decoupled formulation allows exact computation 
of V{t,z), by computation of Vi{t,Xi), value functions of 
lower dimensional spaces. 

Lemma 1: z S J'i(—fo) => z: G Vf G [—T, —to] for 

some to such that 0 < to < T. 

Proof: Suppose z G Fi{—to), then by (|^ we have the 
following: 

• DtV{—to,z) = 0. Thus, V{t,z) becomes independent 
of f at f = —to- 

• since Ti is an open set, there exists a neighborhood 
around z that is contained in .Fi. Thus V{t,z) is also 
independent of f in a neighborhood of z. 

• By ( [TT] ), we have max^^y miud^v D 2 V{—to, z) ■ 
f \z,u,d) > 0. 

Let ti G {to,T] and suppose z £ F 2 at t = —ti. Then, 
by ([^, max„giYmindGD£’zV"(-ti,z) • f{z,u,d) < 0. 
Since / is independent of t, this necessarily means that 
D,V(-ti,z)^D,V{-to,z). 

This implies zJzq in a neighborhood of z such that 
Vi-h,zo)fVi-to,zo). 


V{-ti,zo)-V{-to,zo) = Af0 (12) 

However, by 0, V{t,zo) = V{—to,zo) Vt G (—ti,—fo]- 
In particular, then, we have for any e > 0 

V{-ti,zo) -V{-ti -f e, zq) ^ A 
e e 

This means that VM G K, 3e > 0 such that 

e 

which is a contradiction since V{t,zo) is Lipschitz 
continuous. Therefore, since max^^u laaind^v DzV{t, z) ■ 
f{z,u,d) > 0 Vf < —to, we have that z G Fi{t) Vf G 
[-r,-fo] by ([nj. 


V(t, z) = max Vi(t(z),Xi) Vz G J-'i(f) 

* (15) 

V{t,z) = maxVi{t,Xi) dz G F 2 {t) 

i 

where Vi{t, Xi), i = 1..., N are the viscosity solutions to 
DtVi{t,Xi) + max min D^.V^{t,Xi) ■ fi{x^,Ui,df) = 0 

Ui&Ai diG'Di 

Vi(0,Xj) — 

(16) 

and t{z) is the smallest time such that z G Fi{t), i.e. 

f(z) = inf{z G J'i(t)} (17) 

T>t 

Proof: Case 1: By Lemma [u we have Vz G z G 

J^i(r)Vr < f(z). Therefore, V(t,z) satisfies (|^ Vf < f(z), 
so V{t,z) = V{t{z),z). Case 2 of this proof would then 
imply V{t{z),z) = maxiVi{t{z),x^). 

Case 2: Consider a target set represented by the zero 
sublevel set of the function Z(z), where l{z) = maxili{xi). 
By (|^, we have that 1/(0, z) = max^ li{xi). Define functions 
Vi{t,Xi) such that Vi{0,Xi) = Uixi), then at f = 0, we have 
V{t, z) = maxi Vi{t, xf), and 

DtV{t,z) = ^ l{z = aTgmaxVi{t,x^)}DtVi{t,x^) 

i 

D^V{t,z) = 'V l{i = argmaxI/i(Cxi)}4,L>x.V"i(f,Xi) 

i 

(18) 

where !{•} is the indicator function that is 1 when its 
argument is true and 0 otherwise, and I^. is an matrix in 
of all zeros except for in the rows corresponding to 
the Xi component where it is the identity matrix in 

Now, consider all points z G 7^2(i), in which V{t,z) 
satisfies Substituting V{t, z) = max^ Vi{t, Xi) into 
we have 


EH* = aT:gmaxVi{t,Xi)}[DtViit,Xi)+ 

* ' (19) 

max min Dx,Vi{t,Xi) ■ fi{xi,Ui,di)] =0. 

Equation (|T^ states that in the region where Vi{t,Xi) is 
the maximum among Xj)}^g, we have V{t,z) = 

Vi(t,Xi), where Vi(t,Xi) satisfies |T^ . 

Consider auxiliary functions Wi{t,Xi) which satisfy, for 
all t G [—T, 0] and all Xi, 


Corollary 1: z G F 2 {—t\) => z G 7^2(f) dt G [—fi, 0]. 

Proof: Suppose 3to G [0,fi),z ^ 7^2(—fo) but z G 
d- 2 {—ti)- Since Fi{t) is the complement of F^'it), this 
implies z G 7^i(—<o)- By Lemma ® we must have that 
z G 7^1 (—fi) since ti G [—T, —fo]^ a contradiction. ■ 

We can now state our main theorem. 


DtWi{t,Xi) + max min Dx.Wi{t,Xi) ■ fi{xi,Ui,di) = 0 

Ui &Ai di G'Di 

Wi{0,x^) = li{xi). 

( 20 ) 

By Corollary® we have that Vi(t,Xi) and Wi(t,Xi) both 
satisfy the same PDE with the same terminal conditions, i = 
1... ,N. Therefore, Vi{t, Xi) = Wi[t, Xi),di = 1,..., N. ■ 




C. Decoupled Formulation Algorithm 

Algorithmically, Theorem [T] states the following: 

1) DtV{to,z) = 0 for some to => DtV{t,z) = 0 Vf G 
[-T,to]. 

2) \fzeJ^i{t),DtV{t,z)=0. 

3) \/z G !F 2 {t),V{t,z) = maxj Vi(t, a;^) where Vi{t,Xi) 
satisfies 

This gives us an efficient way to computed V(t, z) by 
computing Vi{t,Xi),i = 1,..., N, transforming the original 
n-dimensional problem of computing V{t, z) into the N rii- 
dimensional problems of computing Xi), i = \ ... ,N. 
Based on the conclusions we drew, the following algorithm 
exactly computes V(t, z), which satisfies (|^, with the above- 
mentioned computation benehts: 

1) Initialize Vi(t, Xi) = li{xi), i = 1 ..., N. 

2) Compute Vi(t, Xi), i = 1..., N, hy solving 

3) Initialize V{t,z) = l{t,z) = maxiVi{t,Xi). 

4) Decrement t from 0 to —T; for each time step t: 

a) Set the auxiliary variable 

V{t, z) ^ max Vi {t, Xi) 

i 

This step correctly computes V(t, z) to be V{t, z) 
for all z G lF 2 (i). 

b) Update the value function 

V{t,z) ^ min{U(f,z),U(f,z)} 

This step correctly computes V{t, z) to satisfy 
DtV(t,z) = 0 for all z G J-'i(f). 


D. Computation Time and Space Complexity Comparison 

For a state space discretization of k grid points in each 
dimension, the computation time complexity decreases from 
0{kF) for the original problem in M", to Oi^^kF') — 
for the N subproblems in ]R"G This is a com¬ 
putation speed improvement of many orders of magnitude. 

Directly solving (|^ on a computational domain S C K" 
has a space complexity of 0{Tk'^), where r is the number 
of time steps of V {t, z) being stored, since we need to store 
an n-dimensional grid for each of the t time steps. The 
algorithm presented in |III-C involves computing Vi{t,Xi) 
on computation domains Xi C = l,...,iV. Each 

Vi{t,Xi) thus has a space complexity of 0{Tk'^'), making 
the overall space complexity "•). 

From Vi{t,Xi), we can then reconstruct U(f, z) in any 
domain Z C A"! x A 2 x... x Xm- Thus, by choosing Z to be a 
small subset of A”! x... x Ajv, we can always avoid additional 
space complexity. This allows us to access V{t,z),z G Z. 
In practice, one would choose Z to be in a small region 
around a state z of interest (eg. the current system state), 
and access the value function V{t, z) as well as its gradient 
Dz(t, z) at z; this allows one to determine of whether z is in 
the reachable set based on the sign of V (t, z), and compute 
the optimal controls u{t),d{t) for both Player 1 and Player 
2 respectively based on DzV{t,z) as follows: 


d* 


argmaxmin£) 2 U(f,z) • f(z,u,d) 

u&A d^X> 

dxg\mnDzV{t,z) ■ f{z,u*,d) 


( 21 ) 


Note that we do not need to store V {t, z) for all z G 
5 = A”! X ... X A'jv- In fact, in many situations, storing 
V{t, z) in the entire S is infeasible since the space com¬ 
plexity is exponential with the dimension of S. With that 
caveat, we restate our algorithm for computing V{t, z) from 
Vi{t,Xi),i = 1,..., N, explicitly noting memory allocation, 
to show that one only needs to store V{t,z) for z G Z C S, 
where Z is a very small subset of S: 


1) Initialize Vi{t,Xi) = li{t,Xi) for Xi G Xi,i = 


2) Compute Vi{t, Xi) in Xi, i = 1..., by solving 16 


3) Initialize V{t,z) = l(t,z) in a small computation 
domain Z C A”! x ... x Xpj. 

4) Decrement t from 0 to —T; for each time step t, 
perform the following computations in Z: 


a) V{t, z) ^ maxi Vi{t, Xi)^ 

b) U(f, z) ^ min{U(f, z),U(f, z)} 


E. Numerical Implementation 

Our proposed decoupled formulation involves solving ( [T6| ) 
for each of the N subsystems. As already mentioned, many 
numerical tools already exist for solving ( [T6| ); we will use 
the implementation in [15]. For the examples in this paper, 
we used the numerical schemes below. 

For the numerical Hamiltonian H{Dx-Vi,Xi) = 
max„g;^ mindgx> DzV{t, z) ■ f{z, u, d) in ( fT^ , we used the 
Lax-Friedrich approximation [26]. For the numerical spatial 
derivatives DanVi{t,Xi), we used a fifth-order accurate 
weighted essentially non-oscillatory scheme [26], [27]. 
For numerical time derivatives DtVi{t,Xi), we used a 
third-order accurate total variation diminishing Runge-Kutta 
scheme [27], [28]. 

Computations were done on a computer with an Intel Core 
i7-2600K CPU running at 3.4 GHz, with 16 GB of memory. 


IV. 4D Quadrotor Collision Avoidance 

Consider a simple quadrotor model consisting of two 
decoupled double-integrators: 


Px — Py — tJy 

Xx - X^XJ Xy - Uy (22) 

U < \Ux\, \Uy\ < U 

Px,Py denote the x- and y-position of the quadrotor, and 
Vx, Xy denote the x- and y-velocity. The control signals 
Ux,Uy are the x- and ^-acceleration of the quadrotor, con¬ 
strained to be between u and u. 

Now, consider two quadrotors in a pursuit-evasion game, 
in which the evader (Player 1) aims to avoid collision, while 
the pursuer (Player 2) aims to cause a collision. The relative 
coordinates of the two quadrotors are given by the following 
state variables: 






5 


V, = [0 1] 


v, = [1 11 


Px,r — Px,i Px,j, Py,r —Pj/,1 Py,2 (23) 

12x,r — '^x,i ^x,j^ ^y,2 

Given the above relative state variables, the relative dy¬ 
namics of the two quadrotors are given by 


Px,r - Py,X - (24) 

'^x^r — ^a:,l ^x,2: ^y^r — ^y,l ^y,2 

Note that this system is decoupled, with xi = 
{Px,r,Vx,r) € as the first decoupled component, and 
X 2 = (Py.r, I'y.r) € as the second decoupled component. 
In the relative coordinates z := {Px,r,Vx^rjPy,r,Vy^r) S M'* 
of the two quadrotors, we define the collision set of size 1, 
representing the configurations in which the two quadrotors 
are considered to have collided, as the following set; 


£ = {z e I < 1} (25) 

with the corresponding implicit surface function l{z) where 
l{z) < 0 ^ z G C. Since we have a decoupled system, let 
Ci,i = 1,2 be the following sets; 


£i = {cci e I Ipx.rl < 1} 

£2 = {X2 G I \Py,r\ < 1} 


with corresponding implicit surface functions li{xi),i = 1, 2. 
Then, we have £ = £in £2 and l{z) = max^ li{xi),i = 1, 2. 

We will set £ as the target set in our reachability problem, 
and compute the backwards reachable set V(f) from £ using 
three methods; 


• Solve directly in to obtain V{t,z), whose zero 
sublevel set represents V{t). 

• Solve ( [Thl ) in K^, i = 1,2, to obtain V{t,z) using our 
proposed decoupled formulation described in Section 

me 

• Compute the analytic boundary of the reachable set. 


For comparison purposes, for the first two methods we 
will compute V{t,z) on the computation domain [—5,5]^. 
However, it is important to recall that for our proposed 
method described in Section III-B| we can significantly 
reduce space complexity by only storing a small part of 
Vit,z). 


A. Reachable Set 

Since the state space of our system is 4D, we visualize 
various {vx,r,Vy^r) slices of the reachable set, whose bound¬ 
ary is given by {z | I4(f = 1.5, z) = 0}. Figures shows 
these slices. The reachable set boundary computed using our 
proposed decoupled method is very close to the reachable 
set boundary computed by solving the full PDF (|^ in 
and to analytic reachable set boundary. Figures zooms in 
on the plots for a closer look. 



Fig. 1; Various {vx,r,Vy,r) slices of the reachable set. 



-1 -0.5 0 0.5 1 -2 -1 0 


v^ = [-1 2] v^ = [-2.5 -1] 



Fig. 2; Various {vx,r, Vyx) slices of the reachable set, zoomed 
in. 


B. Performance 

In order to quantify the computation error, we converted 
V{t,z) into a signed distance function Vsd{t,z) from the 
boundary V{t, z) = 0. This operation was first proposed 
in [29] and can be done by solving the reinitialization 
PDF formulated in [30]; for this operation, we use the 
implementation in [15]. We then evaluated approximately 24 
million analytically-computed reachable set boundary points 
on the Vsd{t, z); the resulting values represent how far each 
of the analytically-computed points are from the numerically- 
computed boundary. The values of Vsd{t,z) on analytic 
boundary points are defined as the computation error. 

Figure shows the error as a function of grid spacing. 
In terms of the maximum error (red curve), the decoupled 
formulation results in a numerically-computed reachable set 
































(vxr, vyr) = (-2, 3) 



Fig. 3: Mean and maximum error of the reachable set 
computed using the decoupled formulation, as a function of 
the grid spacing. 



Fig. 4; Computation time as a function of the number of grid 
points in each dimension. 


boundary that is accurate within the size of the grid spacing 
(black line). On average, the error is approximately an order 
of magnitude smaller than the size of grid spacing (blue 
curve). Furthermore, we can see numerical convergence to 
the analytic solution as the grid spacing size decreases. 

Figure shows the computation time as a function of the 
number of grid points in each dimension. Here, we can see 
that the decoupled formulation is orders of magnitude faster 
than the full formulation, and can be done with many more 
grid points in each dimension. Lastly, the slopes of curves in 
the log-log plot show an 0{k‘^) time complexity for the full 
formulation, and only 0{k^) for the decoupled formulation. 

For the decoupled formulation, when we reconstruct the 
full value function in 4D (blue curve), the computation time 
hardly increases compared to when we do not perform full 
reconstruction (green curve). However, in general, we rec¬ 
ommend that the value function in only a region near a state 
of interest should be computed. Without full reconstruction 
of the value function, we are able to obtain results with many 
more grid points (green curve), improving the accuracy of 
the numerical computation. 

V. 6D Quadrotor Collision Avoidance 

Consider relative dynamics augmented by the velocity of 
evader quadrotor, given in Equation ( |27| ). These dynamics 
are needed to impose a velocity limit on the quadrotor. 

Px,r — '^x,ri Py^T — 

'kx^r — Ha; t^a:; ~ Hy dy (27) 

Ha;,l — Uxt ~ Hy 


(vxr, vyr) = (0, 0) 



(vxr, vyr) = (0, -3) (vxr, vyr) = (2, 2) 



Fig. 5: 2D slices of the 6D reachable set for the augmented 
relative dynamics of two quadrotors. 


For this system, we consider a collision between the two 
quadrotors, as defined previously, to be unsafe configura¬ 
tions. Here, we denote this set Cc 

= {-Z G I \Px,r\ < d, \py^r\ < 2} (28) 

with corresponding implicit surface function lc{z)- 

We also consider configurations in which Player 1 is 
exceeding a velocity limit of 5 in the x- or y- directions 
to be unsafe. This set of configurations is denoted Cs- 

= {-z G I \vx,i I > 5 V |uy,i I > 5} (29) 

We define the target set C to be the union of the above 
configurations; C represents all unsafe configurations. We 
will compute the reachable set V from the target set C as 
follows; 

1) Compute Veit), the reachable set from target set Cq- 

2) Compute Vs{t), the reachable set from target set Cs- 

3) Take the union to obtain V{t) = Vc{t) U Vs(f). 

A. Reachable set 

Taking the union V = Vs{t) U Vc{t), we obtain the 
6D reachable set. To visualize V{t), we compute 2D slices 
of the 6D reachable set at various {vx,r,Vx,i,Vy^r,Vy,i) 
values. This is done without computing the entire 6D reach¬ 
able set by setting the computation domain Z to be in 
a large portion of the {px,nPy,r) plane, at a small range 
of {vx,r,Vx,i,Vy^r,Vy,i) values. The 2D slices at shown in 
Figure]^ Each subplot shows two different pairs {vx,i,Vy^i) 
for a particular {vx,r,Vy^r)- The red boundary represents the 
slice with {vx,i,Vy^i) = (0,0), while the blue boundary 
represents a slice with the evader velocity almost exceeding 
the limit of i; = 5. The blue boundaries contain the red 
ones, because if the evader is already near the velocity limit, 
it would have more limited capability to avoid collisions. 

Figure shows a simulation of the collision avoidance 
maneuver resulting from the reachable set. The red and 
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• -^ Position and velocity 2 
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Fig. 6: A simulation of collision avoidance between the two 
quadrotors using the 6D reachable set. 

blue quadrotors are initially traveling in opposite directions. 
Consider the situation in which the red quadrotor insists on 
staying on its intended path. In this case, the blue quadrotor 
must perform the optimal avoidance control whenever it 
reaches the boundary of the reachable set {t = 1,3), shown 
as the blue dotted boundary. When the blue quadrotor is no 
longer at the boundary of the reachable set (t = 5), it is free 
to perform any control. 

VI. Conclusions and Future Work 

We have presented a decoupled formulation of the HI 
reachability, enabling us to efficiently compute the viscosity 
solution of the HJ PDF (|^, which gives the reachable 
set for optimal control problems and differential games. 
When the system dynamics are decoupled, our decoupled 
formulation allows the exact reconstruction of the solution 
to the full dimensional HJ PDF by solving the HJ PDFs 
corresponding to each of the decoupled components. Our 
novel approach enables the analysis of otherwise intractable 
problems. In addition, our formulation achieves this without 
sacrificing any optimality compared to the full formulation. 
We demonstrated the benefits of our approach using a 4D 
and 6D quadrotor system. 
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